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Abstract 

If the gravitino is the lightest supersymmetric particle and the long-lived next-to-lightest 
sparticle (NSP) is the stau, the charged partner of the tau lepton, it may be metastable 
and form bound states with several nuclei. These bound states may affect the cosmological 
abundances of 6 Li and Li by enhancing nuclear rates that would otherwise be strongly sup- 
pressed. We consider the effects of these enhanced rates on the final abundances produced in 
Big-Bang nucleosynthesis (BBN), including injections of both electromagnetic and hadronic 
energy during and after BBN. We calculate the dominant two- and three-body decays of both 
neutralino and stau NSPs, and model the electromagnetic and hadronic decay products using 
the PYTHIA event generator and a cascade equation. Generically, the introduction of bound 
states drives light element abundances further from their observed values; however, for small 
regions of parameter space bound state effects can bring lithium abundances in particular 
in better accord with observations. We show that in regions where the stau is the NSP with 
a lifetime longer than 10 3 — 10 4 s, the abundances of 6 Li and 7 Li are far in excess of those 
allowed by observations. For shorter lifetimes of order 1000 s, we comment on the possibility 
in minimal supersymmetric and supergravity models that stau decays could reduce the 7 Li 
abundance from standard BBN values while at the same time enhancing the 6 Li abundance. 
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1 Introduction 



The abundances of the light nuclei produced by primordial Big-Bang nucleosynthesis (BBN) 
provide some of the most stringent constraints on the decays of unstable massive particles 
during the early Universe [1-9]. This is because the astrophysical determinations of the 
abundances of deuterium (D) and 4 He agree well with those predicted by homogeneous BBN 
calculations, and also the baryon-to-photon ratio r] = ra 6 / n 7 oc Q b h 2 needed for the success 
of these calculations [10,11] agrees very well with that inferred [12] from observations of the 
power spectrum of fluctuations in the cosmic microwave background (CMB). The value of 
T) = (6.11 ± 0.25) x 10~ 10 that they indicate [13] is now quite precise, reducing one of the 
principal uncertainties in the previous BBN calculations. 

However, it is still difficult to reconcile the BBN predictions for the lithium isotope 
abundances with observational indications on the primordial abundances. The discovery 
of the "Spite" plateau [14], which demonstrates a near-independence of the 7 Li abundance 
from the metallicity in Population-II stars, suggests a primordial abundance in the range 
7 Li/H ~ (1 — 2) x 10~ 10 [15], whereas standard BBN with the CMB value of rj would predict 
7 Li/H ~ 4 x 10~ 10 [10, 11]. In the case of 6 Li, the data [16] lie a factor ~ 1000 above the 
BBN predictions [17], and fail to exhibit the dependence on metallicity expected in models 
based on nucleosynthesis by Galactic cosmic rays [18]. On the other hand, the 6 Li abundance 
may be explained by pre-Galactic Population-Ill stars, without additional over-production 
of 7 Li [19]. 

The concordance between BBN predictions and the observed abundances of D and 4 He 
is relatively fragile and could have been upset by decays of massive unstable particles. Elec- 
tromagnetic and hadronic showers produced in decays occurring during or shortly after BBN 
induce new reactions which may either create or destroy light nuclides. Concrete examples of 
unstable but long-lived particles are found in supersymmetric theories. In particular, models 
based on the constrained version of the minimal i?-parity conserving supersymmetric stan- 
dard model (CMSSM) with a gravitino as the lightest supersymmetric particle (LSP) have 
been considered in this context [20-25]. In these models with a gravitino LSP, the next-to- 
lightest supersymmetric particle (NSP) may be either the lightest neutralino x or the lighter 
stau fi, and the decay lifetime of the NSP can range from seconds to years depending on 
the specific model and gravitino mass. Our previous results concentrated on relatively long 
lifetimes (r > 10 4 s), and the effects of the electromagnetic showers on the light-element 
abundances [4,21-23]. 

The effects of hadronic injections due to late decays of the NSP during BBN have also 
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been studied extensively by other authors [3,6,7,26-28]. In particular, it has been shown 
for relatively short lifetimes of order 10 3 s that decays may simultaneously increase the 
6 Li abundance and decrease the 7 Li abundance [26,27]. On the other hand, it has been 
shown [23] that purely electromagnetic showers cannot reduce the 7 Li abundance sufficiently 
without also overproducing D relative to 3 He [29]. In [20] the authors have calculated the 
principal three-body decays of a stau NSP and, using the results of the analysis [6] have 
explored regions where the 7 Li puzzle can be solved in an unconstrained supersymmetric 
model. The authors of [6] used the PYTHIA [30] model for e + e" annihilation to hadrons 
in order to simulate the hadronic decays of the unstable particle. An analogous simulation 
was used in [24], both to locate regions in the parameter space of the CMSSM which are 
compatible with the BBN constraints, and to solve the lithium problem [25]. 

It has recently been pointed out that, if it has electric charge, the NSP forms bound states 
with several nuclei [31]. Due to the large NSP mass (m nsp ^> m nuc i eo n), the Bohr radii of 
these bound states ~ a ~ lm nucicon ~ 1 fhi are of order the nuclear size. Consequently, nuclear 
reactions with nuclei in bound states are catalyzed, due to partial screening of the Coulomb 
barrier, and due to the opening of virtual photon channels in radiative capture reactions. 
Ref. [31] used analytic approximations to argue that the d(a,7) 6 Li reaction is enhanced by 
an enormous factor ~ 10 6 , leading to 6 Li production far beyond tolerable levels over large 
regions of parameter space. Other effects of bound states were considered in [32,33]. 

Here, we present results from a new BBN code that includes the nuclear reactions induced 
by hadronic and electromagnetic showers generated by late gravitational decays of the NSP, 
together with the familiar network of nuclear reactions used to calculate the primordial 
abundances of the light elements Deuterium (D), 3 He, 4 He and 7 Li. In addition, we include 
the effects of the bound states when the decaying particle is charged. As claimed in [31], these 
bound states lead to large enhancements in otherwise heavily suppressed rates. We calculate 
the abundances of the various bound states and include both the Coulomb enhancement as 
well as virtual photon effects in radiative capture reactions. 

We use as frameworks for this study both the CMSSM and mSUGRA models [22], where 
the NSP could be either the lighter stau or the lightest neutralino. We have calculated the 
dominant two- and three-body gravitational decays of these sparticles and use the PYTHIA 
Monte Carlo event generator and a cascade equation to model the resulting electromagnetic 
(EM) and hadronic (HD) spectra for each point of the parameter space of the supersymmetric 
model. The resulting accurate determinations of the abundances of the light elements, as 
altered by these late injections and the bound-state effects, enable us to delineate regions 
of the parameter space of the supersymmetric models which are compatible with the BBN 
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constraints. In addition, we look for regions where the 6 Li and 7 Li puzzles can be solved in 
the context of these supersymmetric models. We find that for lifetimes r < 10 3 — 10 4 s, the 
enhanced rates of 6 Li and 7 Li production, exclude gravitino dark matter (GDM) with a stau 
NSP. At smaller lifetimes, we see that it is the 7 Li destruction rates which are enhanced, 
facilitating a solution to the Li problems. 

2 Electromagnetic and Hadronic NSP Decays 

The CMSSM is determined by four real parameters, namely the soft supersymmetry-breaking 
scalar mass m , the gaugino mass m^, the trilinear coupling A (each taken to be universal 
at the Grand Unification scale), the ratio of Higgs vevs tan (3, and the sign of the \i parameter. 
Here, for simplicity, we restrict our attention to A = and /i > 0. The mass of the gravitino 
is also a free parameter in the CMSSM, and if is chosen to be less than min(m x ,m f ) the 
resultant model has GDM. In mSUGRA models [22], tan/5 is no longer a free parameter, 
nor is the gravitino mass which is now equal to m at the Grand Unification scale. GDM is 
a inevitable consequence in mSUGRA models when m is relatively small. In this case, the 
NSP is typically the stau, but it may also be the neutralino if A /m < 1.7. The abundances 
of light elements provide some of the most important constraints on this gravitino LSP 
scenario [20-25]. They also impose important constraints on the neutralino LSP scenario, 
since a gravitino NSP would also decay slowly. Here, however, we restrict our attention to 
GDM scenarios with either a stau or neutralino NSP. 

In order to estimate the lifetime of the NSP, as well as the various branching ratios and 
the resulting EM and HD spectra, one must calculate the partial widths of the dominant 
relevant decay channels of the NSP 1 . The decay products that yield EM energy obviously 
include directly-produced photons, and also indirectly-produced photons, charged leptons 
(electrons and muons) which are produced via the secondary decays of gauge and Higgs 
bosons, as well as neutral pions Hadrons (nucleons and mesons such as the K®, ^ 

and 7r ± ) are usually produced through the secondary decays of gauge and Higgs bosons, 
as well (for the mesons) as via the decays of the heavy r lepton. It is important to note 
that mesons decay before interacting with the hadronic background [3,20]. Hence they are 
irrelevant to the BBN processes and to our analysis, except via their decays into photons 
and charged leptons. Therefore, the HD injections on which we focus our attention are those 
that produce nucleons, namely the decays via gauge and Higgs bosons and quark- ant iquark 
pairs. 

1 Analytical results for all of the relevant partial widths will be presented elsewhere [9]. 
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For the neutralino NSP Xi we include the two-body decay channels x ~^ G Hi and 
X —>■ GV, where Hi = h,H,A and V — 7, Z. These are the dominant gravitational decays of 
X, whose analytical expressions have been presented in [21]. In addition, we include here the 
dominant three-body decays x — > G 7* — > G qq, x — > G 7* — > GW + W~, x — > GW + W~ and 
the corresponding interference terms. In general, the two-body channel x G 7 dominates 
the x NSP decays and yields the bulk of the injected EM energy. When the x is heavy 
enough to produce a real Z boson, the next most important channel is x ~~ * G Z, which 
is also the dominant channel for producing HD injections in this case. The Higgs boson 
channels are smaller by a few orders of magnitude, and those to heavy Higgs bosons (H, A) 
in particular become kinematically accessible only for heavy x i n the large-m^ region. 
Turning to the three-body channels, the decay through the virtual photon to a qq pair 
can become comparable to the subdominant channel x ~~ * G Z, injecting nucleons even in 
the kinematical region m x < m 3 / 2 + Mz, where direct on-shell Z-boson production is not 
possible 2 . Finally, we note that the three-body decays to W + W~ pairs and a gravitino are 
usually at least five orders of magnitude smaller. 

For general orientation, we present in Fig. ^ contours (in seconds) for the supersymmetric 
models discussed later and presented in Figs. El and El We display the mi/ 2 ,m planes for 
variants of the CMSSM with different gravitino masses for tan (5 = 10, A = (panels a and 
b), and tan/5 = 57, A = (panel c). Panel d displays an analogous plane for a mSUGRA 
model with Aq/tuo = 3 — a/3 at the GUT scale, in which tan/3 is determined at each point 
by the electroweak vacuum conditions. In addition to the NSP lifetime contours (labelled 
by the lifetime in seconds), we show the boundaries between regions with neutralino and 
stau NSPs (dotted red lines) and the upper limit on the gravitino mass density (solid brown 
lines). These curves are also found in Figs. El and El 

Having calculated the partial decay widths and branching ratios, we employ the PYTHIA 
event generator [30] to model both the EM and the HD decays of the direct products of the 
X decays. We first generate a sufficient number of spectra for the secondary decays of the 
gauge and Higgs bosons and the quark pairs. Then, we perform fits to obtain the relation 
between the energy of the decaying particle and the quantity that characterizes the hadronic 
spectrum, namely dN h /dE h , the number of produced nucleons as a function of the nucleon 
energy. These spectra and the fraction of the energy of the decaying particle that is injected 

2 In principle, one should also include qq pair production through the virtual Z-boson channel \ ~ * G Z* — > 
Gqq [7] and the corresponding interference term. However, this process is suppressed by a factor of Af| 
with respect to \ ~~ * G7* — > Gqq, and the interference term is also suppressed by Mf . Numerically, these 
contributions are unimportant, and therefore we drop these amplitudes in our calculation. 
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m 3 / 2 = 100 GeV , tan (3 = 10 , \i > m 3/2 = 0.2m , tan (3 = 10 , |ll > 




m 1/2 (GeV) m 1/2 (GeV) 

Figure 1: T/ie NSP lifetime contours (in seconds) for the supersymmetric models discussed 
in Figs. [H and\^ 
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as EM energy are then used to calculate the light-element abundances. 

An analogous procedure is followed for the f NSP case. As the lighter stau is predomi- 
nantly right-handed, its interactions with W bosons are very weak (suppressed by powers of 
m T ) and can be ignored. The decay rate for the dominant two-body decay channel, namely 
f -> Gr, has been given in [21]. However, this decay channel does not yield any nucle- 
ons. Therefore, one must calculate some three-body decays of the f to obtain any protons 
or neutrons. The most relevant channels are f — > Gt* — > GZt, f — > Z f* — > GZr, 
f — > tx* — ► GZt and f — » GZt [20]. We calculate these partial widths, and then use 
PYTHIA to obtain the hadronic spectra and the EM energy injected by the secondary Z- 
boson and r-lepton decays. As in the case of the \ NSP, this information is then used for 
the BBN calculation. 

We stress that this procedure is repeated separately for each point in the supersymmetric 
parameter space sampled. That is, given a set of parameters mo, mi/2, A , tan f3, sgn(/x), and 
7713/2, once the sparticle spectrum is determined, all of the relevant branching fractions are 
computed, and the hadronic spectra and the injected EM energy determined case by case. 
For this reason, we do not use a global parameter such as the hadronic branching fraction, 
Bh, often used in the literature. In our analysis, Bh is computed and differs at each point in 
the parameter space. 

3 Electromagnetic and Hadronic Showers During Pri- 
mordial Nucleosynthesis 

The dominant effect of hadronic decays of the NSP during BBN is the addition of new 
interactions between hadronic shower particles and background nuclides 3 . These alter the 
evolution and final values of the light-element abundances, as follows. For each background 
species i, let the rate of interactions of decay hadrons per background particle be I\. Then 
the i abundance per background baryon, Yj = Ui/nB, changes according to 

dtYi = (ftY 4 ) SBBN + (<9^) HAD + (d t Yi) EM , (1) 

where Ysbbn gives the rate of change of the i abundance in standard BBN. We have also 

included in (0) the effects of electromagnetic interactions due to NSP decays, either from the 

decays directly to photons or leptons, or through electromagnetic secondaries in the hadronic 

showers. These are treated as in [4], but are not dominant when hadronic branchings are 

3 The decays of the NSP affect, in principle, the expansion rate of the Universe. However, this effect is 
negligible for NSP abundances low enough to respect the other constraints discussed below. 
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significant. All we need to know is the total EM energy released per decay in any given 
channel, which may become more complicated in the three-body case, but can easily be 
calculated. 

Including hadronic decays in BBN thus amounts to an evaluation of the interaction rates 

(9tYi) HAB = -r^ incl F 4 + Fhb^Yb. (2) 

hb 

The first term accounts for % destruction by hadro-dissociation, where T^nei is the total 
rate (for a given species %) of all inelastic interactions of shower particles with %. The second 
term accounts for production via the hadro-dissociation of heavier background species, e.g., 
Pshowcrttbg — ► d. The sum runs over shower species h and background targets b. In the 
case of the lithium isotopes, production also occurs via the interactions of energetic (i.e., 
nonthermal) mass-3 dissociation products with background 4 He, e.g., 3 Heo; — > 6 ' 7 Li + • • • . 

Consider a nonthermal hadronic projectile species h, with energy spectrum Nh(e,t) = 
dN h /dVde and total number density rih = f N h de. The rate for % production due to hb — > i 
is 

IW-i = J N h (e,t)a hb ^i(e)de. (3) 

The rates thus depend on the decay particle, on the background abundances, and most 
importantly from the point of view of implementation, on the shower development and 
evolution of N h;e in the background environment. 

We wish to follow the evolution of Nh over the multiple shower generations produced by 
the initial hadronic NSP decay products. In the context of BBN, this problem of shower 
development has been approached via direct computation of the multiple generations of 
shower particles [2,3,26]. In this approach, the final particle spectrum is obtained by iterating 
an initial decay spectrum, accounting for both the energy losses and the energy distributions 
of collision products. 

We introduce here an equivalent alternative approach, based on a cascade equation, 
emulating the well-studied treatment of hadronic shower development due to cosmic-ray 
interactions in the atmosphere. The spectrum of h evolves according to 

d t N h {e) = J h (e) - A h (e)N h (e) - d e [b h (e)N h (e)} , (4) 

where Jh is the sum of all source terms, is the sum of all sink terms, and b h is the 
energy-loss rate of particle-conserving processes. The energy-loss term is assumed to be the 
dominant process of energy transfer, in which case tertiary processes are limited to down- 
scattering. The sink term has two contributions, due to elastic and inelastic scattering; 
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the source term has three contributions, due to direct injection, elastic down-scattering and 
inelastic down-scattering. 

Each nonthermal species i evolves according to a cascade equation of the form (J2|, and 
together these constitute a coupled set of equations. Because the source term includes the 
elastic term with iV, inside the integral, these equations are of integro-differential form. The 
integration is therefore not immediate. Previous work on hadronic decays has adopted a 
Monte Carlo approach to the solution; our method is to solve the differential equation. 

We note that because this is an integral equation, we can adopt an iterative approach to 
the solution. To make our initial guess, we ignore the downscattering and solve for Nh with 
decays being the only source. The solution can then be written in terms of the following 
quadrature 



where our initial guess takes = ?x( e ) only. Here the exponential "optical depth" factor 



is a measure of the average number of inelastic interactions over the time taken to lose energy 
electromagnetically from e' to e. 

The full cascades can then be treated iteratively, correcting the approximation to include 
the redistribution and production of nucleons in scattering events. We do this by using 
the previous solution to update the source term J^ +V> by including the downscattering 
terms. These distributions converge after a few iterations. We then insert them into (J2J) and 
solve for the hadro-dissociation rate. This iterative procedure is similar to what previous 
studies have done. However, rather than including the exponent in the integral, they treat 
the exponential R term as a delta function, evaluated at e* (R(e*, e, t) = 1). 

Full details of our method will be given in [9]. 

4 Bound-State Effects 

It has recently been pointed out [31] that the presence of a charged particle, such as the stau, 
during BBN can alter the light-element abundances in a significant way due to the formation 
negatively-charged staus of bound states (BS) with charged nuclei. The binding energies of 
these states are a 2 Zfrrii/2 « 30ZfA{ keV, and the Bohr radii ~ (aZiirii) ~ 1 Zf 1 A^ 1 fm. 
For species such as 4 He, 7 Li and 7 Be, these energy and length scales are close to those of 






(6) 
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nuclear interactions, and it thus turns out that bound state formation results in catalysis of 
nuclear rates via two mechanisms. 

One immediate consequence of the bound states is a reduction of the Coulomb barrier for 
nuclear reactions, due to partial screening by the stau. Since Coulomb repulsion dominates 
the charged-particle rates, all such rates are enhanced. Specifically, for the case of an initial 
state Ai + A 2 , Coulomb effects lead to a exponential suppression via a penetration factor 
which scales as with A = A\Aij{A\ + A 2 ). Introduction of a bound state 

(f, A 2 ) decreases the target charge to Z 2 — 1 and the system's reduced mass number to 
A = A±] both effects lower the Coulomb suppression. We include these effects for all reactions 
with bound states. 

Table 1: Bound-state virtual photon enhancements to radiative capture cross sections for 
various reactions. The third (fourth) column is the threshold energy for the standard (cat- 
alyzed) BBN. The last column is the virtual photon enhancement factor of the catalyzed cross 
section relative to standard BBN. 



EM A(B,j) [X~A} (B, C)X~ Enhancement 

Reaction Transition Qsbbn (MeV) Qcbbn (MeV) ctcbbn/ctsbbn 

d(a,7) 6 Li E2 1.474 1.124 7.0xl0 7 

3 H(«,7) 7 Li El 2.467 2.117 1.0x10 s 

3 He(a,7) 7 Be El 1.587 1.237 2.9xl0 5 

6 Li(p, 7 ) 7 Be El 5.606 4.716 2.9xl0 4 

7 Li(p, 7 ) 8 Be El 17.255 16.325 2.6xl0 3 

7 Be(p, 7 ) 8 B El 0.137 -1.323 N/A 



An additional effect enhances radiative capture channels A 2 (Ai, 7) A by introducing pho- 
tonless final states in which the stau carries off the reaction energy transmitted via virtual 
photon processes. In particular, the 4 He(d, 7) 6 Li reaction, which is suppressed in standard 
BBN, is enhanced by many orders of magnitude by the presence of the bound states, as 
described in [31]. Large enhancements of this type affect other radiative capture reactions, 
notably mass-7 production reactions such as 3 He(a,7) 7 Be and destruction reactions such as 
7 Li(a, 7) n B. We have included these as well: the corresponding enhancement factors appear 
in Tabled 

Bound-state formation and reaction catalysis occurs late in BBN. The binding energy 
E hin for the [f, 4 He] bound states is 311 keV, for [f, 7 Li] 952 keV, and for the [f, 7 Be] 1490 
keV. The latter are quite high, of order nuclear binding energies, and indeed the large 7 Be 
binding plays an important role in forbidding 7 Be destruction channels that otherwise would 
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be energetically allowed. The capture processes that form these bound states typically 
become effective for temperatures T c « E^- m /30] this means that 7 Be states form prior to 7 Li 
states, with 4 He states forming last. At these low temperatures one can ignore the standard 
BBN fusion processes that involve these elements. 

To account for bound state effects, an accurate calculation of their abundance is necessary. 
To do this we solve numerically the Boltzmann equations (13) and (14) from [32], that control 
these abundances. If X denotes the light element, and ignoring the fusion contribution as 
described before, the system of the two differential equations for the light-element and bound 
state abundances can be cast into the form 

Yx = { j^r{Y x n,-Y BS n' 1 ) 

Y BS = -Y x , (7) 

where Y B s,x — nBs,x/s and n? is the stau number density. The thermally-averaged capture 
cross section (a c v) and the photon density n'^ for E > E^, ia , are given in Eqs (9) and (15) 
in [32], respectively. H is the Hubble expansion and dot denotes derivatives with respect to 
the temperature. As initial condition, we assume that the bound state abundance is negligible 
for a temperature of a few times T c . In our numerical analysis we solve the system (j7J) for 
X = 4 He, 7 Li, 7 Be to obtain the corresponding Ybs at temperatures below T c . We assume 
that the bound state is destroyed in the reaction. That is, we do not include additional 
bound-state effects on the final-state nuclei such as 6 Li. 

As we see in the following section, bound state effects indeed greatly enhance 6 Li produc- 
tion as found in the analysis of d(a,7) 6 Li by [31]. Our systematic inclusion of bound state 
effects finds that 7 Li is also significantly altered. The most important rates are for radiative 
capture reactions, which enjoy large boosts due to virtual photon effects. In particular, bound 
state 7 Li production is dominated by the 3 H(a,7) 7 Li and 3 He(o;,7) 7 Be rates. Destruction is 
dominated by the channel with the lowest Coulomb barrier, namely 7 Li(p, 7)2 4 He. Note that 
7 Be destruction channels are less important, since mass-7 is largely in 7 Li at T > 60 keV, 
and because the high binding energy of [f, 7 Be] makes [f , 7 Be] + p — > 8 B + f energetically 
forbidden with Q = -1.3 MeV (see Table H}. 

Interestingly, the bound state perturbations lead to net 7 Li production in some parts of 
the parameter space of the models we study, and net destruction in others. Net production 
occurs when the stau is sufficiently long-lived (r f > 10 4 sec) that 4 He bound states are abun- 
dant enough to drive bound state enhanced production stronger than bound state enhanced 
7 Li destruction. On the other hand, within a window of slightly shorter lifetimes, staus per- 
sist long enough to form 7 Li bound states, but then decay before forming 4 He bound states. 
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This leads to net 7 Li destruction. Thus we see that 7 Li is quite sensitive to the r properties, 
and potentially offers a strong probe of the existence and nature of bound states. 

5 Results and Discussion 

As described earlier, we work in the context of the CMSSM or mSUGRA. Our primary 
goal is to examine the effect of bound-state interactions on the final abundances of the 
light elements. To this end, we display a selection of results for specific supersymmetric 
planes both with and without the effect of bound state interactions. All results shown fully 
incorporate the effects of electromagnetic and hadronic showers. A more complete selection 
of results as well as constraints which go beyond the MSSM will be presented in [9]. 

We begin with results based on CMSSM models with A = 0, fi > and tan (3 = 10. We 
display our results in the (m^jfno) plane, showing explicit element abundance contours. In 
Fig. EK, we show the element abundances that result when the gravitino mass is held fixed 
at 777.3/2 = 100 GeV in the absence of stau bound state effects. To the left of the near- vertical 
solid black line at 7771/2 — 250 GeV, the gravitino is the not the LSP, and we do not consider 
this region here. Immediately to the right of this line is a red dot-dashed line. To the left 
of this, the Higgs mass is below the current experimental bound 114 GeV. The diagonal red 
dotted line corresponds to the boundary between a neutralino and stau NSP. Above the line, 
the neutralino is the NSP, and below it, the NSP is the stau. Very close to this boundary, 
there is a diagonal brown solid line. Above this line, the relic density of gravitinos from NSP 
decay is too high, i.e., 

^^tt NS ph 2 > 0.12. (8) 

777AT5P 

Thus we should restrict our attention to the area below this line. Note that we display the 
extensions of contours which originate below the line into the overdense region, but we do 
not display contours that reside solely in the upper plane. 

We start with the solid orange line labelled 3 He/D = 1. To the left of this curve, the 
3 He/D ratio is greater than 1, which is excluded [29]. For small mo, this excludes gaugino 
masses less than about 1100 GeV, which is similar to the result found in [23]. To the right 
of this curve, the ratio of 3 He to D is acceptable. The very thick green line labelled 7 Li = 
4.3 corresponds to the contour where 7 Li/H = 4.3 x 10~ 10 , a value very close to the standard 
BBN result for 7 Li/H. It forms a 'V shape, whose right edge runs along the neutralino- 
stau NSP border before shooting up at 777 1 / 2 ~ 3900 GeV. Below the V, the abundance of 
7 Li is smaller than the standard BBN result. However, for relatively small values of 777 1 / 2 , 
the 7 Li abundance does not differ very much from this standard BBN result: it is only 
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Figure 2: Some (m^jfflo) planes for A = 0, jj, > and tan/3 = 10. In the upper (lower) 
panels we use m 3 / 2 = 100 GeV (777,3/2 = 0.2 m ). In the right panels the effects of the stau 
bound states have been included, while in those on the left we include only the effect of the 
NSP decays. The regions to the left of the solid black lines are not considered, since there 
the gravitino is not the LSP. In the orange (light) shaded regions, the differences between 
the calculated and observed light- element abundances are no greater than in standard BBN 
without late particle decays. In the pink (dark) shaded region in panel d, the abundances lie 
within the ranges favoured by observation, as described in the text. The significances of the 
other lines and contours are explained in the text. 
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when mi/2 ?t 3000 GeV that 7 Li begins to drop significantly. This is seen by the additional 
(unlabeled) thin green contours showing 7 Li/H = 3 x 10~ 10 (solid), and 2 x 10~ 10 (dashed). As 
can be seen in Fig. [TJi, the the stau lifetime drops with increasing mj/a, and when r ~ 1000 
s, at rax 12 ~ 4000 GeV, the 7 Li abundance has been reduced to an observation-friendly value 
close to 2 x 10~ 10 as claimed in [27]. 

However, for this case with m 3 / 2 = 100 GeV the 6 Li abundance is never sufficiently high 
to match the observed 6 Li plateau for the same parameter values where 7 Li is reduced. The 
6 Li/ 7 Li ratio is shown by the solid blue contour labeled 6 Li/ 7 Li = 0.15. Note that there is 
also a small contour loop at this value at small mo centred around m^j ~ 1000 GeV. Inside 
the loop, the lithium isotope ratio is acceptable, which is also the case to the right of the 
nearly vertical contour at large mx/%. At large mj/2, the contour for 6 Li/ 7 Li = 0.01 is shown 
by the thin blue line. To the right of this contour, including the region where 7 Li ~ 2 x 10~ 10 , 
the 6 Li abundance is too small. 

Finally, we show the contours for D/H = 2.2 and 4.0 xlO -5 by the solid purple contours 
as labeled. The D/H = 2.2 xlO -5 contour is a small loop within the 6 Li/ 7 Li loop. Inside 
this loop D/H is too small. Between the two curves labeled 4.0, the D/H ratio is high, but 
not necessarily excessively so. 

In summary, the acceptable regions found in Fig. EH break down into 2 areas: one between 
the two loops labeled 2.2 and 0.15 and to the right of the 3 He/D = 1 line, where D/H is 
larger than 2.2 xl0~ 5 and 6 Li/ 7 Li < 0.15. However, in this region, the 7 Li abundance is 
very similar to the standard BBN result, which may be considered too high. Alternatively, 
one could consider very large m\ji where once again D/H < 4.0 x 10~ 5 . Here, 7 Li is in fact 
acceptably low, but the 6 Li abundance is far below the plateau value. As a better illustration 
of our results, we have shaded these two regions. The orange (lighter) shaded region is where 
3 He/.D < 1, 6 Li/ 7 Li < 0.15, 2.2 x 10~ 5 < D/R < 4.0 x 10" 5 and 7 Li/H < 4.3 x 10~ 10 . 

Turning now to Fig. Et>, we show the analogous results when the bound-state effects are 
included in the calculation. The abundance contours are identical to those in Fig. EH above 
the diagonal dotted line, where the NSP is a neutralino and bound states do not form. We 
also note that the bound state effects on D and 3 He are quite minimal, so that these element 
abundances are very similar to those in Fig. EH- However, comparing panels a and b, one sees 
dramatic bound-state effects on the lithium abundances. The loop of 6 Li/ 7 Li = 0.15 centred 
about m 1( / 2 = 1000 GeV has now gone due to the large abundance of 6 Li produced by bound- 
state catalysis. Indeed, everywhere to the left of the solid blue line labeled 0.15 is excluded. 
In the stau NSP region, this means that > 3000 — 3500 GeV. Moreover, in the stau 
region to the right of the 6 Li/ 7 Li = 0.15 contour, the 7 Li abundance drops below 9 x 10~ n 
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(as shown by the thin green dotted curve) and D/H > 4 x 10 -5 for m.1/2 < 3500 — 4000 GeV. 
Only when m\i% > 3500 — 4000 GeV does the D/H abundance drop back to acceptable levels 
with good abundances for 7 Li, but 6 Li is now too small to account for the plateau. Thus, 
for a constant value of m 3 / 2 = 100 GeV, the bound-state effects force one to extremely large 
values of mi/ 2 primarily due to the enhanced production of 6 Li, as shown by the orange 
shaded region. For this value of the gravitino mass, there are no regions where both lithium 
abundances match their plateau values. 

We do not display the results for m 3 / 2 = 10 GeV, but the bound-state effects (and the 
results) are less dramatic. Without the bound-state effects included, the 6 Li abundance is 
generally too small, while the 7 Li abundance is very similar to standard BBN in the stau 
NSP region. The gravitino relic density is a factor of 10 smaller in this case and some of 
the neutralino NSP region is allowable. In the neutralino NSP region, D/H is too high 
unless m.1/2 ^ 2000 GeV. At (m.1/2, mo) — (2100,1000), there is a region where D/H and 
7 Li are acceptably small, though 6 Li/ Li is very small. The bound-state effects again set a 
lower limit on m 1( / 2 in the stau NSP region in this case. When m.1/2 ^ 1300 GeV, both Li 
abundances drop and approach their standard BBN values. Once again, in no region are 
both lithium isotopes at their plateau values. 

It is also interesting to consider cases in which the gravitino mass is proportional to m . 
In Fig. |2fc, we fix m 3 / 2 = 0.2 m and neglect the bound-state effects. The choices of contours 
are similar to those in panels a and b. The gravitino relic density constraint now cuts out 
some of the stau NSP region at large m.1/2 and large mo, but allows a small neutralino NSP 
region at low m.1/2. As before, we are constrained to the right of the curve labeled 3 He/D = 
1, though in this case the constraint is not very strong in the stau NSP region. The 7 Li/H = 
4.3 x 10~ 10 contour again forms a 'V shape and one is restricted to lie below the 'V. In most 
of the stau NSP region, 7 Li remains relatively high but begins to drop at large ~ 3000 
GeV as r approaches 0(1000) s (see Fig.Q). The region where the 6 Li/ 7 Li ratio lies between 
0.01 and 0.15 now forms a band which moves from lower left to upper right. Thus, as one 
can see in the orange shading, there is a large region where the lithium isotopic ratio can 
be made acceptable. However, if we restrict to D/H < 4.0 x 10~ 5 , we see that this ratio is 
interesting only when 7 Li is at or slightly below the standard BBN result. However, we do 
note that as one approaches the gravitino density limit at (mi/ 2 ,m ) — (4100,1000), it is 
possible to have 6 Li/ 7 Li ~ 0.06 and 7 Li/H ~ 2.3 x 10" 10 at the expense of D/H ~ 6.2 x 10~ 10 . 

The bound-state effects when m 3 / 2 = 0.2 mo are shown in Fig. EJi. Once again, we see 
that the increased production of both 6 Li and 7 Li excludes a portion of the stau NSP region 
where mi/ 2 ^ 1500 GeV for small mo- The lower bound on mi/ 2 increases with mo- In this 
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case, not only do the bound-state effects increase the 7 Li abundance when m.1/2 is small (i.e., 
at relatively long stau lifetimes), but they also decrease the 7 Li abundance when the lifetime 
of the stau is about 1500 s. Thus, at (mi/ 2 ,m ) ~ (3200,400), we find that 6 Li/ 7 Li ~ 0.04, 
7 Li/H ~ 1.2 x 10" 10 , and D/H ~ 3.8 x 10 5 . Indeed, when m 1/2 is between 3000-4000 GeV, 
the bound state effects cut the 7 Li abundance roughly in half. In the darker (pink) region 
(which has no analogue in the other panels), the lithium abundances match the observational 
plateau values, with the properties 6 Li/ 7 Li > 0.01 and 0.9 x 10~ 10 < 7 Li/H < 2.0 x 10~ 10 . 

For a larger ratio of m 3 / 2 = thq, the gravitino relic density forces us to relatively low 
values of m < 500 GeV. As in the case described above, the viable region at low mi/ 2 
is excluded by the bound-state effects, and we find increased 7 Li production due to bound 
states at high mm. Qualitatively, this case is similar to that when m 3 / 2 = 0.2 mo, though 
most features are compressed to lower values of m . 

In Fig. El we show some examples of results from CMSSM models with tan/? = 57, 
and mSUGRA models. The dominant effect of increasing tan (3 is on the neutralino and 
stau relic densities. At low tan /3, the relic density of the neutralino is generally high except 
along a narrow strip where neutralino-stau co- annihilations are important and yield a density 
with in the WMAP range. At large tan/3, new annihilation channels are available. Most 
predominant is the the s-channel annihilation of neutralinos through the heavy Higgs scalar 
and pseudoscalar, causing large variations in the relic density across the plane, particularly 
at large 7711/2 and tuq. These variations have an impact in GDM scenarios, as the abundance 
of decaying particles varies. 

In Fig. EK, we show the (mu2, m ) plane for tan /3 — 57 (which is near the maximal value 
for which the electroweak symmetry breaking conditions can be satisfied), and m 3 / 2 = 0.2 mo- 
The dark green shaded region at very low mi/ 2 is excluded by b — > 37 decays. Notice that 
the constraint on the gravitino relic density (shown by the solid brown line) no longer tracks 
the neutralino-stau NSP border. At m^j ~ 1200 GeV, it shoots upwards towards large m . 
This is due to the s-channel annihilation pole (where 2m x = hia) which decreases the relic 
density. Consider now the behavior of the 7 Li abundance as mi/ 2 is increased at a fixed value 
of m = 2000 GeV. At small m.1/2 ^ 1000 GeV, the neutralino is the LSP and results are not 
shown. When tan/5 = 10, the relic neutralino density is very large, and when mi/ 2 ^ 2400 
GeV the lifetime is greater than 5 x 10 4 s, and 7 Li destruction process are very efficient. At 
larger m^, the lifetime decreases and hadronic production effects begin to dominate, and 
the 7 Li abundance becomes very large. When tan/5 = 57, the relic density of neutralinos 
is 2-3 orders of magnitude smaller when m-i/ 2 ^ 1500 GeV, for the same value of mo- As a 
result, 7 Li destruction is suppressed. As mi/ 2 is increased, and we move away from the pole, 
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Figure 3: Some more (mi/2,mo) planes for \x > 0. In the upper panels we use m 3 / 2 = 0.2 mo 
and tan/3 = 57, whilst in the lower panels we assume mSUGRA with m 3 / 2 = m and 
Aq/ttiq = 3 — v^3 as in the simplest Polonyi superpotential. In the right panels the effects of 
the stau bound states have been included, while in those on the left we include only the effects 
of the NSP decays. As in Fig. [U the region above the solid black line is excluded, since 
there the gravitino is not the LSP. In the orange shaded regions, the differences between 
the calculated and observed light-element abundances are no greater than in standard BBN 
without late particle decays. The meanings of the other lines and contours are explained in 
the text. 
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the neutralino density increases, dropping the 7 Li abundance for long lifetimes. As for lower 
tan/?, as m^j is further increased and the r x decreases, the 7 Li abundance becomes large 
again, until one hits the neutralino-stau NSP border. The 'V'-shaped 7 Li contour at large 
mi/2 is visible here as well, though it appears squeezed as the NSP border is moved up at 
large tan (3. 

Just below the NSP border, we see another distinctive feature in Fig. EH- The 6 Li/ 7 Li 
ratio, which is generally too large when the neutralino is the NSP, drops dramatically inside 
a narrow diagonal strip. This occurs because the annihilations of staus are here dominated 
by a similar s-channel pole. Inside this strip, the density of staus is very small, and element 
abundances approach their standard BBN values. At lower mo, over much of the plane with 
a stau NSP, the 7 Li and D abundances are close to their standard BBN values, while 6 Li is 
enhanced. In this case, there is a substantial orange shaded region where the light-element 
abundances are no less acceptable than in standard BBN. 

Our results for tan/5 = 57 and m 3 / 2 = 0.2 m when the bound-state effects are included 
are shown in Fig. EJd. For tan f3 = 57, the stau lifetimes are somewhat longer than the 
corresponding lifetimes when tan/5 = 10. This means that the bound-state effects are 
apparent over a larger portion of the plane with a stau NSP. Both lithium isotope abundances 
are significantly higher. Without the bound states, the 7 Li abundance varies little from 
its standard BBN value, but with their inclusion the 7 Li abundance is somewhat higher, 
generally about 5 x 10~ 10 . The effect on 6 Li is larger. Without the bound-state effects, the 
6 Li/ Li ratio remains small unless either TO1/2 and/or itlq are relatively large. Even then, the 
ratio only increases to a few percent unless toi/2 ~ 3500 — 5000 GeV and m ~ 1200 — 1900 
GeV, where it exceeds 0.15. With the inclusion of the bound states, in much of the stau NSP 
region the 6 Li abundance is too high, exceptions being the area where s-channel annihilation 
occurs or in the lower right corner of the displayed plane. There is no region where the 
light-element abundances lie in the favoured plateau ranges. 

Finally, we come to an example of a mSUGRA model. Here, because of a relation between 
the bilinear and trilinear supersymmetry breaking terms: Bq = Ao — mo, tan/3 is no longer a 
free parameter of the theory, but instead must be calculated at each point of the parameter 
space. Here, we choose an example based on the Polonyi model for which A /m = 3 — \/3. 
In addition, we have the condition that = mo- In Fig. Ob, we show the mSUGRA model 
without the bound states. In the upper part of the plane, we do not have GDM. We see 
that 3 He/D eliminates all but a triangular area which extends up to m = 1000 GeV, when 
m i/2 — 5000 GeV. Below the 3 He/D = 1 contour, D and 7 Li are close to their standard BBN 
values, and there is a substantial orange shaded region. We note that 6 Li is interestingly 
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high, between 0.01 and 0.15 in much of this region. 

As seen in Fig. Eli, when bound-state effects are included in this mSUGRA model, both 
lithium isotope abundances are too large except in the extreme lower right corner, where there 
is a small region shaded orange. However, there is no region where the lithium abundances 
fall within the favoured plateau ranges. 

6 Conclusions 

We have calculated in this paper the cosmological light-element abundances in the presence 
of the electromagnetic and hadronic showers due to late decays of the NSP in the context of 
the CMSSM and mSUGRA models, incorporating the effects of the bound states that would 
form between a metastable stau NSP and the light nuclei. Late decays of the neutralino NSP 
constrain significantly the neutralino region, since in general they yield large light-element 
abundances. The bound-state effects are significant in the stau NSP region, where excessive 
6 Li and 7 Li abundances exclude regions where the stau lifetime is longer than 10 3 — 10 4 s. 
For lifetimes shorter than 1000 s, there is a possibility that the stau decays can reduce the 
7 Li abundance from the standard BBN value, while at the same time enhancing the 6 Li 
abundance. A more complete account of our calculations will be given in [9], where more 
examples of CMSSM and mSUGRA parameter planes will be presented, and the possibility 
of matching the favoured lithium abundances will be discussed in more detail. 
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